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Abstract 

Physical concepts developed to describe instabilities in traffic flows 
can be generalized in a way that allows one to understand the well- 
known instability of supply chains (the so-called "bullwhip effect"). 
That is, small variations in the consumption rate can cause large vari- 
ations in the production rate of companies generating the requested 
product. Interestingly, the resulting oscillations have characteristic fre- 
quencies which are considerably lower than the variations in the con- 
sumption rate. This suggests that instabilities of supply chains may 
be the reason for the existence of business cycles. At the same time, 
we establish some link to queuing theory and between micro- and 
macroeconomics . 

1 Introduction 

Concepts from statistical physics and non-linear dynamics have been very 
successful in discovering and explaining dynamical phenomena in traffic flows 
PQ. Many of these phenomena are based on mechanisms such as delayed 
adaptation to changing conditions and the competition for limited resources, 
which are relevant for other systems as well. This includes pattern formation 
such as segregation in driven granular media [2] and lane formation in colloid 
physics [3 J or biological physics (pedestrians, ants) Other examples are 
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clogging phenomena at bottlenecks in freeway traffic [3] , panicking pedestrian 
crowds [HI, or granular media [7j. In the following study, we will focus on the 
phenomenon of stop-and-go traffic jH] and its analogies. 

Recently, economists and traffic scientists have wondered, whether traffic 
dynamics has also implications for the stability and management of supply 
chains [HI EH] or for the dynamics of business cycles jTTj. To explain busi- 
ness cycles, many theoretical concepts have been suggested over the decades, 
such as the Schumpeter clock ^2]. These are usually based on macroeco- 
nomic variables such as investment, income, consumption, public expendi- 
ture, or the employment rate, and their interactions. In contrast, Witt et al. 
|TT] have recently suggested to interpret business cycles as self-organization 
phenomenon due to a linear instability of production dynamics related to 
stop-and-go waves in traffic or driven many-particle systems. In order to il- 
lustrate their idea, they have transferred a continuous macroscopic traffic 
model and re-interpreted the single terms and variables. 

The author believes that this is a very promising approach to understand 
business cycles, but instead of simply transferring macroscopic traffic models, 
suggests to derive equations for business cycles from first principles, which 
means to derive the dynamics on the macroscopic level from microscopic 
interactions in production systems. This would also make some contribution 
to the goal of understanding macroeconomics based on microeconomics (or, 
in a wider sense, based on the "elementary interactions" of individuals, here: 
production managers). 

In order to make some progress in this direction, we will generalize [12] 
an idea suggested by Daganzo to describe the dynamics of supply chains 
|n]. Like the work by Armbruster et al. [TP] , his approach is related to traf- 
fic models as well, but he focusses on models in discrete space in order to 
reflect the discreteness of successive production steps. In order to describe 
the non-linear dynamics of production processes or interrelated economic 
sectors, we will have to generalize these ideas to complex supply networks. 
In Sec. El we will first discuss "macroscopic" business cycles in a sector ally 
structured economy and compare them with stop-and-go traffic. Afterwards, 
in Sec. HH we will develop a more fine-grained, "microscopic" description of 
the management of dynamically interacting production units and relate it to 
classical queueing theory (which mainly focusses on stochastic fluctuations 
of production processes in a stationary state). Later on, in Sec. H] we will 
construct a mathematical relation between the microscopic and the macro- 
scopic level of description, while some further research directions and other 
potential applications are indicated in the outlook of Sec. |3 
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2 Modelling "Macroscopic" Supply Net- 
works 

2.1 Economic Production Sectors 

We will first investigate a simplified economic system with U production 
sectors B generating certain kinds of products I G {1, . . . , P}. The market 
for products J as a function of time t shall be represented by the stock 
level ("inventory") Nj(t). Let us assume that, in each production cycle, the 
production sector B generates p B products of kind / and requires c J B products 
of kind J ("educts") for production. If Qsif) is the number of production 
cycles per unit time, i.e. a measure of the production rate, productivity or 
"throughput" of sector B, the quantity of products generated per unit time 
is p B QB(t), while the quantity of educts consumed per unit time is c J B Q B {t) 
(see Fig.Q). The temporal change of the quantity of products / in the market 
is, therefore, given by the conservation equation 

^ = £(Ps - c B )Q B (t) . (1) 

The production rate Qs(t) will be specified later on in Sec. El For the time 
being, we will assume Qb{1) to agree with the actual feeding ("arrival") rate 
\ B discussed in the next paragraph: 

\ B (t) = Q B {t). (2) 



2.2 The Feeding Rates 

In the absence of capacity constraints, we would have the relation X B = Pb^b 
for the time-dependent feeding rate Ab, where /is is the processing rate of 
sector B, i.e. the potential production rate in the case of no inefficiencies. 
p° B reflects the desired utilization. A value of p° B = 0.7 is reasonable, as it 
guarantees a relatively high production rate at moderate and reliable waiting 
times for finishing the products (see Sec. IH.lj) . However, the actual utilization 
Pb = ^b/i^b does not necessarily agree with the desired one, p° B . In order to 
illustrate this, let us assume that T B is the average time needed to deliver to 
production sector B products of kind J, which are required for production 
(which will be called "educts" in the following). If Nj is the quantity of 
actually available products of kind J, the resulting maximum quantity of 
products that can be delivered per unit time is Nj/T B . As c J B educts are 
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Figure 1: Schematic illustration of the virtual and actual product flows for 
the case of two kinds of products, I and J. The picture displays the various 
model variables related to production sector B, which is represented by the 
oval structure in the center. Rectangles correspond to markets, arrows to 
product flows. 

required per production cycle, one cannot have a higher production rate 
than Nj/(T B c B ) per unit time, which defines the maximum feeding rate \ J B 
of educt J into production unit B: 

\ J B = Nj/(T J B c J B ) . (3) 

Moreover, as the products can only be finished, if all required educts are 
available in a production cycle, the actual feeding rate is given by the 
minimum of these values and of the desired feeding rate p° B fiB (see Fig. [TJ): 

X B (t) =mm(p Bf x B ,{\ J B }). (4) 

If c J B = 0, we set — * oo, so that the corresponding term does not have any 
impact on the value of the minimum function (fl|). 

Note that production has some analogy with chemical reactions, where 
one also requires a certain quantity of educts to produce other (chemical) 
products. For chemical reactions in threedimensional space, however, the 
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reaction rate is given by multiplication of (a power of) the chemical concen- 
trations of the educts, i.e. the above minimum function is replaced by an 
algebraic product. This difference is comparable to the Probabilistic AND 
and the Fuzzy AND in Fuzzy Logic, which are applied when several condi- 
tions are to be met at the same time. The Probabilistic AND is based on 
a multiplication of the logical values, while the Fuzzy AND corresponds to 
their minimum. 



2.3 Adaptation of Production Speeds and Transport 
Capacities 

One important aspect is the adaptation of production to the time-dependent 
inventories Nj(t). The adaptation of the desired utilization p° B is delayed (see 
Sec. EJ), and a change of the processing rate fis, which requires an adaptation 
of production sectors or capacities, takes a considerable time as well. We will 
assume 

^ = i[W B ,...)-W (5) 
at tb 

with a typical adaptation time tb- Wb{^-/Zb, • • • ) reflects the production rate 
of sector B in steady state as a function of the inventories Nj. For the time 
being, we will assume 

1 _ STPb*! 



Z B (t) ^ Nj(t) 

i.e. the production speed is adapted only to the inventories Nj of the products 
J generated by sector B, for which p B > 0. X B is an additional and constant 
prefactor. This formula is defined in a way which gives always positive values 
Zb($), if not all coefficients p J B are zero and the Nj(t) stay positive (see 
below). Normally, the production speed Wb{1/Zb, ■ ■ ■ ) will increase with 
decreasing inventories Nj, but it saturates due to financial or technological 
limitations and inefficiencies. In this study, we will assume 

WsWZ* ...) = max (a b 1 + ^Jff^ , o) , (7) 

where Zb could be called "scaled inventories" and As, Bb, and Db are 
suitable parameters (see Fig. El). 

Furthermore, we will assume here that the transport capacities are 
adapted in parallel with the production speed and represent the propor- 
tionality constants by V B . This implies 

Ab = P>BViN.j/c J B , (8) 
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Figure 2: Steady-state production speed Wb as a function of the scaled in- 
ventories Zb for the parameter Bb = 0.2, arbitrary Ab > 0, and various 
values of Db- The production speed decreases with increasing inventories. 
Throughout this paper we use the parameter values Bb = 0.2 and Db = 8. 



corresponding to the specification 1/Tg = P° b PbV^, which results in the 
following formula for the feeding rates: 

X B (t) = p B // B min fl, {^^}) • (9) 
Herein, the minimum extends over all indices J. 



2.4 Bull- Whip Effect and Stop-and-Go Traffic 

For a linear supply chain with c B = 8b,i+i and p B = Sb,i (where Sk,l = 1, if 
K = L, and 5k, l = otherwise) we obtain the particular set of equations 

^ = Qiit) - Q I+1 (t) = \j(t) - A J+1 (t) . (10) 

In the capacity-constrained case (characterized by small values of Vg), this 
leads to 

^ = V^N^t) - V I+1 (t)Ni(t) (11) 
Vj(t) = p^jVj- 1 (12) 

= ± [Wjil/Zr) - fa] (13) 
6 



with 
and 



or 

^ = 1 [y/^i/Zj) - v#)] . (14) 

Interestingly, Eqs. ([TTj) and (|14j) basically agree with the macroscopic traffic 
flow model by Hilliges and Weidlich 14J, where (jllj) is analogous to the 
equation for the vehicle density and (fT4"j) corresponds to the equation for 
the average vehicle speed Vj in street segment I. These equations behave 
linearly unstable with respect to perturbations of Nj(t), if 77 exceeds a certain 
threshold (see Fig. OJ), which depends on the maximum slope W'il/Zj) of 
W(l/Zj) [Hj. As drivers (over-)react with a time delay to a changing traffic 
situation in front, stop-and-go traffic can emerge. The frequently observed 
instability of supply chains, called the "bull-whip effect", occurs for similar 
reasons (e.g. in the "beer distribution game" (XSJUEl), see Fig. El 




Figure 3: Minimum and maximum inventories as a function of the adaptation 
time 7~b. The other model parameters are specified as in Fig. 03 but the pertur- 
bation amplitude has been chosen five times smaller. The difference between 
both curves is the amplitude of the bullwhip effect, i.e. the time-dependent 
variation in the inventories. Note that there is a critical adaptation time, 
below which perturbations are not amplified. In this case, the investigated 
linear supply chain behaves stable. 

The mechanism behind this instability is the delay 77 in the adaptation 
of the production speeds and transport capacities, which implies an over- or 
under-production. The repeatedly or periodically resulting high inventories 
are due to temporary bottlenecks in the supply chain and could be avoided by 
appropriate control functions W B = W B {{Nj(t)}, {dNj/dt(t)}, . . . ) fT7|ITH]. 
Note that, instead of unstable production with high inventories and low pro- 
duction speeds, one may reach the same average, but stable throughput at 
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low inventories and high production speeds (see Fig. Hj). However, this is nor- 
mally related with higher energy and maintainance costs, so that production 
tends to operate in the linearly unstable regime. 
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Figure 4: Steady-state throughput as a function of the scaled inventories 
Zb for the parameter Bb = 0.2, arbitrary Ab, and various values of Db 
in the capacity-restricted case of low transport rates V^. For large enough 
values of Db, the curve has a maximum at finite inventories. In this case, high 
steady-state flows Qb appear twice: for low (and linearly stable) inventories 
and higher (but potentially unstable) inventories (see also Fig. EJ). 

If the transport capacities are not a limiting factor (i.e. the parameters 
Vb are large), instead of fTTj) we have the set of equations 

^ = P?(^/-P? + i(*W- (15) 

Together with this basically corresponds to a particular microscopic 

traffic model, the so-called optimal velocity model [ID] , if we restrict our 
comparison to the linear regime around the steady state and identify p®fii 
with the actual velocity of vehicle J, but Nj with the inverse distance to the 
next car ahead (apart from proportionality constants). It is known fTf\ ITH] 
that this model is linearly unstable, if 

V I I - 1 W' I (1/Z I )>^-, (16) 

i.e. if the adaptation time 77 exceeds a certain threshold which depends on 
the slope W&I/Zt) of Wiil/Zj). 
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2.5 Business Cycles 

The most interesting point is the reaction of the system to a perturbation in 
the throughputs for example a periodic perturbation of the so-called 

consumption rate Qu+i(t) (see Sees. 14. HI and l4.4j) . The resulting oscillations 
in the inventories can be much slower and are usually synchronized among 
different production sectors (see Figs. EU and fT2"j) . Therefore, they may 
explain business cycles as a self-organized phenomenon with slow dynamics 
on the time scale of several years. 

250 i ■ 1 1 1 1 




1000 2000 3000 4000 5000 
Time t 

Figure 5: Variations of the inventories Nj(t) of U = P = 10 economic 
sectors /, triggered by a constant consumption rate with a slight peri- 
odic perturbation, namely Vu(t) = V A 1 1 1 °W / /(1)[1 + 0.1 sin(O.lt)]. One can 
clearly see the emergence of slow and synchronized oscillations of the in- 
ventories Ni(t), which are triggered by small and fast perturbations in the 
consumption rate Vn(t). The unit time is one day and the model param- 
eters are A B = Aj = 100 /V, Vf' 1 = V = 10" 4 , and r B = tj = 90. 
The initial and boundary conditions are iVj(0) = No(t) = 20 = Xj, 
V}(0) = VWi(X}/Nj(0)) = VWj(l), and V n (t) (see above). 

The reduction in the resulting oscillation frequency compared with the 
perturbation frequency is also known from stop-and-go traffic. It has been 
explained by the non-linearity in the model equations. However, there is a sig- 
nificant difference between the dynamics of traffic flows and supply chains: 
While stop-and-go waves have a characteristic amplitude independently of 
the average vehicle density, the amplitude of oscillations in the inventories 
change continuously in the capacity-constrained case (see Fig. In other 
words, the phase transition from stable to unstable traffic flow is hysteretic 
(i.e., of first order) while the phase transition from stable to unsta- 
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ble supply chains in the capacity-constrained case appears to be continuous 
(i.e., of second order). This seems to be a particular property of the Hilliges- 
Weidlich model, while the optimal- velocity model mentioned in Sec. 12.41 is 
known to display a hysteretic transition ^Hj- Therefore, a hysteretic transi- 
tion is found for supply chains which are not contrained by their transport 
capacities. The transition between these two different regimes would be in- 
teresting to investigate. 



o 



> 

3 

| 
3 





S 

3 






Figure 6: Minimum and maximum inventories as a function of the scaled 
steady-state inventory Zj(0). The other model parameters are specified as 
in Fig. but the perturbation amplitude has been chosen five times smaller 
(of the order of 2%). The difference between both curves reflects the am- 
plitude of the time- dependent variation in the inventories. Note that, in the 
capacity-constrained Hilliges-Weidlich case investigated here, the change of 
the amplitude (and the associated kind of non-equilibrium phase transition) 
is continuous. The supply chain is stable with respect to perturbations, where 
both curves agree. This is the case at small and large steady-state inventories. 

Note that, in order to find the emergence of slow oscillations, i.e. of busi- 
ness cycles, we do not need to have a linear supply chain. Supply networks 
can display similar features (see Sec.|3J). The only requirement is that no sta- 
tionary state exists or the stationary state is linearly unstable with respect 
to perturbations. This is intensively studied for particular network types in 
a recent scientific collaboration [18J. 

Summarizing our present insights, preconditions for the emergence of 
business cycles are 

• large values of the adaptation times tb, 
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• significant changes of the production speeds Wb with changing inven- 
tories, 

• the presence of perturbations, if a steady system state exists, 

• no consideration of forecasts or poor forecasts of the future time devel- 
opment of the inventories Nj(t), 

• a highly non-linear interaction among the production sectors. 

We should finally underline that the above considerations can be applied to 
economic systems with any number of sectors and any network structure. 
We do not even need to assume a sectorally structured economy, as the same 
kinds of equations apply on the "microscopic" level of production networks, 
only with a significantly higher number of equations (see Sec. l4.1|) . The special 
cases discussed above have been chosen only for illustrative reasons. 

3 "Microscopic" Model of Production Pro- 
cesses 

In this section, we will formulate a generalized model for the dynamical in- 
teraction of production processes. Compared to the description of interacting 
economic sectors in the previous section, this approach may be called micro- 
scopic. Note that the concept developed in this section is required for two 
reasons: First in order to be able to describe real production processes, which 
involves buffers and other variables and calls for a more complex model. Sec- 
ond in order to allow the derivation of the macroeconomic dynamics from 
microeconomic assumptions regarding the production management of single 
companies (see Sec. HJ). Readers not interested in these aspects may skip this 
section and continue with Sec. |SJ 

3.1 The Production Units in Terms of Queueing The- 
oretical Quantities 

We will now investigate a system with u production units (machines or facto- 
ries) b G {1, 2, . . . , u} producing p products i, j G {1, 2, . . . ,p}. The respective 
production process is characterized by parameters c£ and p\: In each produc- 
tion cycle, production unit b requires c£ products ("educts") j G {1, . . . ,p} 
and produces p\ products i G {1, . . . ,p}. The number of production cycles 
of production unit b per unit time is a measure of the throughput and shall 
be represented by Qb{t)- For production in the steady state, we can often 
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approximate this quantity in terms of variables from queueing theory 
For this, let X b be the feeding ("arrival") rate, C b the number of parallel pro- 
duction channels, //& the overall processing ("departure") rate (i.e. C b times 
the processing rate of a single production channel), 

p b = \ b /nb (17) 

the utilization, and S b the storage capacity of production unit b. Then, the 
following relationships apply (see Fig. Ej) : 

Q b {t) = (1 - p x b )\ b {t) = (1 - rf)n b {t) . (18) 

The functions p b and p b reflect a reduction in the efficiency of the feeding 
(arrival) and the processing (departure) in production unit b. They depend 
on p b , Cb, Sb, and possibly other quantities as well. In some cases, p b is the 
probability of rejecting arrivals, when the storage capacity Sb is reached, i.e. 

pt(pb,Cb,S b ) = Pb(Sb), (19) 

where P b (l) is the probability of having a queue of length l b in production 
unit b. On the other hand, the average number of active production channels 
is 



E 

h— h—Gi 



kPbih) + °b p M • (20) 



Consequently, the relative reduction p b in the processing rate due to {C b — l b ) 
empty channels is 

P^( Pb , C b , S b ) = ^ P b {h) • (21) 

In more complex production systems, the factors p b and p^ would be gener- 
alized measures of inefficiencies in the production process, which would be 
functions of all relevant process parameters. 

According to Little's Law Lb = XbTb for queueing systems, where A& = 
(1 — Pb)K is the average arrival rate, we can express the throughput also in 
terms of the average queue length L b (pb, Cb, Sb) and the average waiting time 
Tb(pb, Cb, Sb): 

= L b (p b , C b , S b ) , 22 , 
Tb(pb, Cb, Sb) 

For example, for a system with infinite storage capacity S b — > oo, we have 

Q b = min(A b , /i b ) = fi b min(p b , 1) . (23) 
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Figure 7: Schematic illustration of a production unit b as a queueing system 
with a limited storage capacity S b and Cb parallel production channels. The 
arrival rate A&, the departure rate pb, and effects of inefficiencies are indicated. 



For a M/M/l : (S^/FIFO) process (one channel with first-in-first-out serv- 
ing, storage capacity Sb, Poisson-distributed arrival times and exponentially 
distributed service intervals), one finds for \b < fib 

Note that both, the expected value and the standard deviation of the queue 
length and the waiting time diverge for p b — > 1. Therefore, efficient produc- 
tion is normally related to pb < 0.7 [22] • 



3.2 The Feeding Rates of Production Units 

In the absence of capacity constraints, we just have the relation X b = p b p b 
for the feeding rate, where the actual utilization p b agrees with the desired 
utilization p b , e.g. p° = 0.7. However, a lack of required educts may lead to 
a reduction of A;,. In that case, the feeding rate is limited by the minimum 
arrival rate A\ of required educts j, divided by the quantity of educts 
needed for one production cycle. We will assume that the maximum arrival 
rates A b are given by p b p b V b P b) where P b denotes the quantity of educts 
stored in the input buffer, p b p b V b is the maximum transport rate due to 
capacity contraints V b for getting educt j from the input buffer into the 
production unit b, and the prefactor p b p b suggests that the transport capacity 
is adapted to the production speed (generalizations are possible). Therefore, 
if A^(t) = PbPf)V^J^(t)/c|, denotes the maximum feeding rate for educt j into 
production unit b, the actually resulting feeding rate is 

A fe (t)=min(pO/i fe ,{A J J), (25) 

where the minimum extends over all indices j. The product flows are illus- 
trated in Fig. |H1 
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Figure 8: Schematic illustration of the product flows together with the model 
variables related to production unit b, which is represented by the oval struc- 
ture in the center. Rectangles correspond to buffers, arrows to product flows. 

When we assume that the flows S 3 ab of products j to b from the delivering 
production units a can be fed into the production process in parallel (in- 
stead of having to go through the input buffer first), we have the generalized 
relationship 



(26) 



with e = 1. The previously discussed case (requiring delivery through the 
input buffer) corresponds to e = 0. In any case, the fraction of the supply 
^2 a S J ab , which is not needed to satisfy the production requirements c^Qb, is 
delivered to the input buffer. 



3.3 Input and Output Buffers 

Let us assume each production unit b has input buffers for required prod- 
ucts i ("educts") and output buffers (a "warehouse") for the products. We 
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will assume the input buffers are filled with Il(t) educts i G {1, . . . ,p} and 
the output buffers with 0\{t) products waiting to be delivered. If S l ab {t) de- 
notes the delivery flow ( "supply" ) of products % from production unit a to b, 
the change of an input buffer stock with time is given by the conservation 
equation 

^ = £*&(*) (27) 

a 

as J2 a Sab(t) i s the quantity of products i delivered from various sources (pro- 
duction units) a, and c\Qb is the quantity of educts i used up for production 
per unit time (see Fig. |HJ). Analogously, the dynamics of an output buffer 
stock is determined by the equation 

^=p'Mt)-J2 S L(t), (28) 

c 

as p{Qb is the quantity of newly generated products j, and ^ c S bc (t) are 
deliveries to other production units c (see Fig.|HJ). The following specifications 
in this paper ensure the non- negativity conditions I b (t) > and 0\(t) > 
(see also Sees. IH.5I and l4.1|l . We will, for example, assume that the delivery 
flows S l ab are basically given by the order flows ("demands") D l ab , but limited 
by the maximum transport rates V* b p b [ib for delivering the available products 
from output buffer O l a to production unit b: 

Si b (t) = minfc(t), Ol(t)V:M . (29) 

If required, suitable specifications can also guarantee I b (t) < J*' max and 
O l b {t) < (9* >max with maximum input and output buffer stocks /*' max and 
0^' max . We may, for example, multiply S l ab by [1 — (^/^' max ) K ], which stops 
the delivery flow of product i when the respective input buffer of production 
unit b is full. Analogously, we may multiply p°/i& by Ylji^ ~ (^b/^&' max ) K ]> 
which stops the production when one of the output buffers of b is full. High 
values of n produce a hard cutoff, while small values of k > 1 can describe 
cases, where the production efficiency goes down even before the buffer size 
is fully used up. 

3.4 Adaptation of Production Speeds and Transport 
Capacities 

One important aspect is the adaptation of production to changing demand. 
On the one hand, a change of the processing rate fib is expensive and time 
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consuming. On the other hand, the adaptation of the desired utilization p° 
is delayed by the average queuing time T b and other factors. We will assume 

^^ = ^[W b (l/Z b ,...)-^] (30) 

with a typical adaptation time r& significantly greater than the adaptation 
times of the other variables. (That is, why we do not consider possible time 
delays in the other mathematical relations, here.) Wb(l/Zb, . . . ) is some con- 
trol function reflecting the strategy by the production manager in adapting 
the utilization and/or processing rate to the output buffer stocks J a or other 
variables. In the following, we will assume 



(31) 



i.e. the management strategy is only sensitive to the perceived stock levels N b 
of the products j generated by unit b, for which p ] b > 0. X 3 b are proportionality 
constants. The perceived stock levels will be specified by 

N J b =0i + 5j20i. (32) 

For 5 = 0, the management takes into account only the own output buffer 
stock 3 h , while for 5 = 1, it tracks also the output buffer stocks 0{ of com- 
peting production units a. Normally, the control function Wb will increase 
with decreasing perceived stock levels N b , but the function Wb(l/Zb, . . . ) sat- 
urates due to financial, spatial, or technological limitations and inefficiencies 
in the processing of high order flows. Again, we will use a function of the 
form 

«'-^-(*ITwtw' G ) (33) 

with suitable parameters A b , B b , and Db- Ab corresponds to the maximum 
production speed. 



3.5 Order Flows, Delivery Networks, and Price Dy- 
namics 

The flow of orders is basically given by the flow c\Qb of educts required for 
the production by unit b. Deviations can be reflected by a correction function 
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Ub(Il, 11°), where If denotes the desired input buffer stock. If the proportions 
q l ab with ^ a q l ab = 1 reflect orders from different producers a, we have 

D l ah = qIAQMlin ) (34) 

and, therefore, 

S l ab (t) = mm[qi b ctQ b U b (li if ), (35) 
(see Eq. (129)) ). When the first term is the smaller one, this implies with (127)1 

dP 

= J2^4Q b u b (H,n°)-4Q b 

a 

= [U b (P b , if ) - MQ b {t) . (36) 

The function U b {P b , If) should be chosen greater than 1, if P b is smaller than 
the desired input buffer stock If, smaller than 1 for P b > If, and otherwise 
1, e.g. 

U b {Il,lf) = lf/Il (37) 

Without this correction function, the input buffer tends to be emptied in the 
course of time. 

Together with c\ and p\, the fractions q l ab characterize the delivery or sup- 
ply network. One can imagine various scenarios. For example, the fractions 
could be assumed constant or modeled by an evolutionary selection equation 
with a selection rate v |2~3*] : 



dt V aby "' ^ 



FL(t)-YFUt)qUt))qUt)- (3* 



The specification of the fitness depends on the relevant parameters. For ex- 
ample, it could be treated constant. However, in some cases, it makes sense 
to relate the fitness F % ab to the inverse of the real or virtual costs ("prices") 
p l ab of product i, when delivered from production unit a to b: 

Fi b {t) = l/pl b (t) . (39) 

Assuming a law of supply and demand, one conceivable specification would 
be 

where p° are the costs when the supply S l ab agrees with the demand D l ab . For 
example, a reduction in the price level would require a slight generalization 
of Eq. flU. 

Other specifications are, of course, possible as well, as the above formulas 
partly depend on the strategies of the human decision makers involved. 
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3.6 Calculation of the Cycle Times 

Apart from the productivity or throughput Q b of a production unit, pro- 
duction managers are also highly interested in the cycle time, i.e. the time 
interval between the beginning of the generation of a product and its comple- 
tion. Let us first discuss the process cycle time t b between entering the queue 
of production unit b and leaving it, assuming that all c\ required educts i 
for one production cycle are transported together and located at the same 
place in the queue. The problem is similar to determining the travel times of 
vehicles entering a traffic jam. 

According to Eq. (|22|). the average waiting time is determined as the 
quotient T b = L b /Q b of the average queue length L b and the throughput Q b , 
if production operates in the steady state. Let us now generalize this formula 
to situations in which the inflow 

Q? = (l-p*)\ b (41) 

into production unit b is time-dependent and possibly differs from the time- 
dependent outflow 

Q? t = (l-rf)l* (42) 

(see Fig. EJ). In reality, this time-dependence results from fluctuations in the 
production process and breakdowns of machines, etc. In some cases, one can 
use the length-dependent formulas 

#<« = { others <«) 

and 

«-(« = { £fc/cJt<S w 

(see Sec. 13. We will now derive a delay-differential equation for the cycle 
time under varying production conditions. For this, let l b {t) be the actual 
length of the queue in production unit b at time t. The change of this length 
in time is given by the difference between the inflow and the outflow at time 
t: 

^ = Q?(t)-Qr(t)- (45) 

If the production unit b has C b channels, an educt entering the production 
queue at time t must move forward l b (t) — C b steps, before it is finally served 
by one of the C b channels. If, after entering the queue at time t, tl(t) denotes 
the waiting time until one of the channels is reached, and if the average 
processing rate of a single channel is /J, b /C b , the average serving or treatment 
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time is given by C b / fib(t + tl). The overall time £& required for the processing 
of the product is, therefore, given by the sum of the waiting time ££ and the 
treatment time by one of the channels: 

*<« = « () + ^fm (46) 

(which replaces the average value T b ). On the other hand, the waiting educts 
move forward Q b nt steps per unit time. For this reason, the waiting time ££ (£) 
is given by the implicit equation 

t+t* b {t) t+t* b {t) t 

hit) -C b = J dt' Q™\t') = j dt 1 Q™\t') - J dt' Q™\t') . (47) 

t — oo — oo 

Identifying the time-derivative of this equation with Eq. ()45|) results in 



QT(t) - Qr(t) = = QT(t + t* b (t)) \i + ) - QT lt (t) , ( is) 

which leads to the delay-differential equation 



dt Qr(t+n(t)) l - (49) 

As the production initially starts with a waiting time of i£(0) = (when 
the factory or production unit h is opened), this equation can be solved 
numerically as a function of the outflow Q1 ut (t'). In this way, it is possible 
to determine the waiting time t* b (t) and process cycle time t b . In the future, 
approximation methods shall be developed for cases where ££(0) is not known 
or Q b nt (t') cannot be controlled or anticipated. A rough approximation would 
be to replace these values by mean values (cf. Eq. 

In a similar way, one can calculate the waiting time t\ in the input buffer 
I b , resulting in 

< = Z a s ab (t) _ 

dt dQr(t+ti(t)) • 1 ; 

The waiting time t\ in the output buffer 3 h is given by 

dt J2 c sL(t + ti(t)) 

Finally, the transport time £^ fe between the output buffer of production unit 
a and the input buffer of production unit b can be estimated by 

4 = ^ ■ (52) 

V ab 
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The total cycle time is, then, calculated as the sum of the respective waiting, 
serving, and transport times. 

4 Relation between Production Processes 
and Macroeconomics 

We will now give some support for the "macroscopic" model of supply net- 
works applied in Sec. El by relating it to the "microscopic" dynamical pro- 
duction model developed in Sec. El 



4.1 Just-in-Time Production 

Let us start with deriving a simplified model of production processes. For 
this, we will summarize the input and output buffers, defining the stock 
levels ( "inventories" ) of markets for the products % by 

b b 

According to the balance equations (J27|) and (J2%|) . we find 

^ = £(pj - 4)Q b (t) . (54) 

b 

Moreover, typical for just-in-time production, we will assume negligible input 
buffers 

Jj = ij° = 0, which implies iV t = ^Oj. (55) 

b 

Consequently, for e = 1 and V* b = V b , the feeding (arrival) rates become 

\ b (t)= P y b mm[l,{^^} ) . (56) 

Furthermore, for 5 = 1 we obtain N b = ] a = Nj and the related adap- 
tation equations 



dt r b 



, my 



(57) 
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We will now assume Qb = A&, which is usually given for small enough 
channel utilizations pb or sufficient storage capacities S&. In that case, one 
would have to solve Eq. (f57jl together with 

Note that formula ()58)1 maintains non- negative inventories Ni(t), as required. 
To show this, let us assume that t would be the first point in time where some 
inventory Nj vanishes, say Ni(t) = 0. The inventory Ni(t) could only become 
negative for dNi/dt < 0, which would require c l b > 0. But then, the minimum 
in formula (|5H)) would be V b l Ni{t) / c\ = 0, which contradicts our assumption 
and proves non- negativity (at least for V b 3 > 0). 



4.2 Micro-Macro Link 

In the following, we will try to reduce the number of equations by an aggrega- 
tion procedure, which keeps the structure of the above model equations. Let 
us first define production sectors B by summarizing those production units b, 
which are characterized by the same adaptation times Tb and a proportional 
throughput Qbif). With suitable constants tb and kb, we can then assume 

n = t b and Qbif) = k b QEi{t) ■ (59) 

The total throughput, processing and feeding rates of production sector B 
are defined by 

QB = S ^Qb, HB = ^2fi>b, A B = ^A 6 , (60) 

beB beB beB 

where b € B indicates that b belongs to production sector B. The equations 
for the inventories keep their structure 

^ = E^-4)QB(t), (61) 



dt u 



if we define 



P i B = J^Pb k b alld C B = Yj^ kb ( 62 ) 

beB beB 

with kb = Qb(0)/Qs(0) and use the relation J2b = J2s^2beB- Moreover, in 
order to get the equations 

X B (t) = p B p B mmUI^^^\] (63) 
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and 



dt tb 



MY, 



p j b x b 



PbPb 



we need to set 



Vi = Vi<J B /4 and X B =Xlp>/p 



B 



(64) 



(65) 



for all j (requiring that the resulting values on the right-hand sides depend 
only on B, but not on b G B), and 



p°b = E p>»/»b , w B (i/z B ) = ]T w b (i/z b 



(66) 



&eB 



with i/z B = EjAxi/Nj = ^jrixi/Nj = 1/Z b . 

We may also summarize all those markets i which show a proportional 
dynamics of Ni(t) in time. This assumes 



Ni(t) = XiNtf) , 
if we denote the proportionality factors by Xi and define 



(67) 



(6£ 



where j 6 / indicates that % is part of the market sector I. Introducing 



Pb = ^2pb and c B = ^ 

i6/ ie7 



(69) 



and summing up Eq. (|61J) over i G / yields the further reduced set of equa- 
tions 

(/A/ E^-4)g B (t). (to) 

B 



dt 



In order to obtain the equations 



Ai?(i) = p^Bmin f 1 ) | 



'B 



and 



d{p%PB) 1 



Pb/^B 



(71) 



(72) 
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we have to define 



and 

X " ^ £i*-V (?4) 

Based on assumptions (}59|) . (}65j) . and (jBTj) . one could achieve a consider- 
able reduction in the number of equations, leading from a microscopic model 
of production processes involving management decisions to a macroscopic 
model of interacting economic sectors. Conditions (|5^j). (JfiHj) . and (|6T|) pre- 
suppose similar production processes of the summarized production units 
b and proportional coefficients regarding the summarized markets i. While 
our aggregation method can be generalized to time-dependent parameters 
kb{t) and Xi(t), we still require u = tb, and the resulting values of V^c^/c^, 
X^p J b /p B should only depend on B, but not b G B. These are the main criteria 
for defining production and market sectors based on empirical data. 



4.3 Dynamic Input-Output Model 

In Eq. (I70j) . the sum over B extends from 1 to U + 1, with B — U + 1 
representing the final consumer sector. Splitting it up, the resulting equations 
read 

= E(Pfl - 4)Qb(*) - <& +1 Qtf + i(f) , (75) 

B=l 

as = 0. Similarly, J in formula (J5fij) runs from to P with J = 

representing a market sector for basic resources. Therefore, 

M0=AM,*-(l,^4^})><>, (76) 

where J e {1, . . . , P}. 

Let us now concentrate again on the case Qs(t) = A#(t) and define the 
consumption rate from market sector / by 

Yjit) = <i +1 Q u+1 (t) . (77) 

When we define the U production sectors B through the respective kind 
of products they produce, we can set p B = ksSsj (see Fig. EJ). Without 
restriction of generality, we can choose fee = 1, which just defines the unit 
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quantity in which we measure products of market sector /. In the stationary 
case with dNj/dt = for all /, we then obtain the equation 

u u 

c 'bQb + Yi = J2 PbQb = Qi , (78) 

B=l B=l 

which corresponds to Leontief's input-output model of macroeconomics [21]. 
Therefore, the above model of supply networks can be considered as dynam- 
ical generalization of this classical economic model. 

We should note that, in view of the instability of supply networks, i.e. 
the existence of the bullwhip effect and of business cycles, the applicability 
of steady-state concepts in economics is questionable. People have, therefore, 
tried to formulate dynamical input-output models for a long time in order to 
take into account investment strategies and other aspects. However, many of 
these approaches have turned out to be inconsistent or useless (cf. Ref. |23]). 
In contrast, the above dynamical input-output model results naturally from 
a much more general model of supply and production networks. Investment 
strategies could, for example, be taken into account by a suitable specification 
of Wb{- ■ ■), which may not only be chosen as a function of the inventories 
Nj, but also of time derivatives such as dQi/dt, dYj/dt, or dNj/dt. 



4.4 Specification of the Boundary Conditions 

In the following, we will discuss some examples. If we assume something like 
a conservation of materials or value, this implies certain constraints, which 
reduce the potential complexity of the related supply networks. First of all, 
the quantity of products / consumed by the production sectors B should be 
generated somewhere, i.e. 

J2 c b = J2pb- ( 79 ) 

B B 

Second, the quantity of educts consumed by some production sector B cor- 
responds to the quantity of its generated products, i.e. 

E c * = E^- ( 8 °) 

For closed systems (with c° B = and c{j +l = 0), this implies that the sum 
over all inventories Nj is constant, i.e. 

E^r = E^(t)E(^-4) = o. (si) 
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Figure 9: Top: Schematic illustration of the actual product flows according 
the the input-output model with p B = 5b,i- The circle summarizes the pro- 
duction sector (oval) and its market (rectangle). The associated product flow 
is represented by a thick arrow, the product flows into the production sector 
by thin arrows. Bottom: Particular case of a linear supply chain, in which 
each production unit receives goods only from the previous one. 

In the following, we will again discuss the case p B = 5b,i- In this particular 
case, relations (|8())l and (f?H|) imply the normalization conditions 

^4 = 1 = ^4. (82) 

B I 

As in Sec. I4.3| we will extend the sums over I from to U and the sums over 
B from 1 to U + 1. In this way, we define 

4 = 1-^4 and 4 +1 = l-]T4. (83) 

I B 

Values 4 > allow us to describe the inflow of basic resources 1 = 0, while 
4+i > allows one to describe the depletion of products by a consumer 
sector B = U + 1, the aging of products, or the generation of products of 
minor quality. The boundary conditions are completely defined by specifying 
No(t) and Qu+i(f), see Eqs. fTSj) and (|75)l. This treatment is fully consistent 
with the intuitive one of linear supply chains in Sec. 12.41 
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5 Impact of the Supply Network's Topology 



We will now discuss simulations based on the equations specified in Sees. 14.11 
to 14.41 for the three different supply networks sketched in Figs. ITUk -c. each 
with five levels: (a) a linear supply chain with 5 production units, (b) a 
"supply ladder" with 10 production units, and (c) a hierarchical supply tree 
with 31 production units. 

By introducting random variables which were assumed to be 
equally distributed in the interval [ — 77, 77] , we have taken into account 
a heterogeneity rj in the individual parameters characterizing the dif- 
ferent production units. Here, we have chosen N (t) = N = 20 = 
X£, JVj(O) = 20(1 + £/), r B = 180(1 + £b), V° = V = 0.05, 
Vi = Vc J B (l + C J b) for J > 0, fP B (0)nB = W B (20/N B (0)), Q u+1 (t) = 
W u+1 {20/N ) min(l, AT z (t)y/4 +1 )[l + 0.1 sin(0.04t)] and A B = 100/V. 

For the linear supply chain, we have c B = 1, if I delivers to b, otherwise 
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For the supply ladder and the hierarchical supply tree, we have c B = 0.5, if 
an arrow points from i to b (see Figs. ITUb . c), otherwise c B = 0. Specifically, 
for the supply ladder we have 
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Figure 10: Illustration of different supply networks: (a) linear supply chain 
with five production units, (b) "supply ladder" with five levels, and (c) hier- 
archical supply tree. Circles summarize production sectors and their markets, 
as in Fig. El while arrows represent the product flows among the sectors. 
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Figure 11: Time-dependent inventories for a supply network with five layers 
and identical model parameters (rj = 0). The dynamics is the same for all 
three network topologies displayed in Fig. i-e. for the linear supply chain, 
the supply ladder, and the hierarchical supply network. 



while for the supply hierarchy, we have 
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c° B and c{j +1 are defined in accordance with Eq. (jHSJ)- These specifications 
guarantee that, for 77 = 0, i.e. if the production units are characterized by 
identical parameters, the dynamics of the inventories is the same for all 
three discussed network topologies (see Fig. 111)). However, the topology mat- 
ters a lot, if we have a heterogeneity rj > in the model parameters (see 
Figs. I12fe-f). As our dynamical model of supply networks assumes non-linear 
interactions, changes of rj can have large effects. The same applies to small 
changes in N (see Fig. [HJ) or in the relaxation times tb (see Fig. |3J). 
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Figure 12: Time- dependent inventories for supply networks in the case of 
heterogeneous parameters (with a heterogeneity of r\ = 0.2). The dynamics 
drastically depends on the topology of the supply network: (a) linear supply 
chain, (b) supply ladder, (c) hierarchical supply network. Compared with 
Fig. one can conclude that heterogeneity in supply networks can consid- 
erably decrease the undesired oscillation amplitudes in the inventories. The 
strongest effect by far is found for supply ladders, which is relevant for the 
design of robust supply networks. 
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6 Summary and Outlook 



In this contribution, we have sketched a dynamical theory of supply net- 
works. Here, we could only present first results and indicate possible future 
research directions, which may contribute to the interdisciplinary field of 
econophysics [213 • The proposed theory is developed to help understand the 
complex non-linear phenomena in production and supply networks, in partic- 
ular their breakdowns, instabilities and inefficiencies. We have also sketched 
how to derive "macroscopic" equations for the dynamics of a sectorally struc- 
tured economics from "microscopic" equations describing single production 
processes, involving strategical decisions of production managers reflected by 
the control functions W^. The resulting model is a dynamical generalization 
of the classical input-output model of macroeconomics. For the particular 
case of a linear supply chain, one can relate it with the Hilliges-Weidlich 
model, which has originally been developed for traffic flow. These equations 
can describe the "bullwhip effect" due to their linear instability in a certain 
regime of operation. The underlying mechanism is the slow adaptation of 
the processing rate to changes in the order flows or stock levels in the mar- 
ket. Interestingly, the resulting oscillations in the inventories of the different 
products % have a characteristic frequency, which can be much lower than 
the underlying fast variations in the consumption rate. Depending on the 
network structure, these oscillations synchronize among different economic 
sectors and may explain business cycles as a self-organized phenomenon with 
slow dynamics. These conclusions are not restricted to linear supply chains, 
but can be generalized to many other supply networks, which are linearly un- 
stable with respect to perturbations. In reality, business cycles are, of course, 
less regular and of smaller amplitudes than in Figs. [3] or ^2 However, the 
above model allows to reflect these aspects in a natural way by inclusion of 
fluctuations, heterogeneity, additional capacity constraints, and realistic net- 
work structures. Investigations with empirical data and for particular kinds 
of networks are on the way. It should also be noted that already deterministic 
models of supply networks can show irregular, non-periodic behavior such as 
chaotic dynamics [131211, which calls for suitable control concepts [2Ej. 

Compared to traffic dynamics, economic and production systems have 
some interesting new features: Instead of a continuous space, we have dis- 
crete production units, and the production speed as a function of the inven- 
tories is different from the empirical velocity-density relation in traffic. Due 
to the minimum condition (|76|h production systems may operate in differ- 
ent regimes, and small changes of parameters may have tremendous effects. 
For example, we may have a transition from small oscillations of relatively 
high frequency to large oscillations of low frequency. Apart from this, the 
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management strategies can vary to a large extent, and with this the control 
functions . . ). With suitable strategies, the oscillations can be mitigated 
or even suppressed (HllIZlllHl- Moreover, production systems are frequently 
supply networks with complex topologies rather than linear supply chains, 
i.e. they have additional features compared to (more or less) one-dimensional 
freeway traffic. They are more comparable to street networks of cities. 

Apart from some equations which were not further applied in this study, 
most of the proposed model equations were conservation equations, equations 
given by the product flows, or relations derived with stochastic concepts used 
in queuing theory. They reflect the transport and interaction of products, so 
that the physics of driven many-particle systems and of complex systems 
can make some significant contributions to the new multidisciplinary field of 
self-organization phenomena in production and supply networks. 

Future work will have to address questions such as the relevance of the 
network structure for the resulting dynamics, possible control strategies, the 
role of the market and pricing mechanism, etc. This particularly concerns the 
specification of the function Wf,, and the equations suggested in the paragraph 
on delivery networks and price dynamics. Here, we have already seen that the 
supply network's topology and the level of heterogeneity in production sys- 
tems have a significant impact on the resulting dynamics: Parameter changes 
can have tremendous effects. For example, a heterogeneity in the parameters 
characterizing the different production units can stabilize the production and 
the market considerably, while in traffic flow, heterogeneity has normally a 
destabilization effect. The stabilization of supply networks through hetero- 
geneity probably could explain why the variations in economic systems ap- 
pear to be less dramatic than in our simulations, but additional inefficiencies 
and capacity restrictions (corresponding to finite values of Sb, /*' max ; and 
0^' max ) are probably another reason. This and the role of heterogeneity for 
the micro-macro link will be studied in more detail in the future. Should it 
turn out that the micro-macro link requires a high degree of homogeneity in 
the parameters, it would be favourable to simulate economic dynamics based 
on a microscopic model of production and supply networks in the future. 
The tendency for globalization certainly increases the degree of homogene- 
ity, but it also tends to generate larger oscillations in the inventories, i.e. 
more serious over- and underproduction. At least in some markets, there are 
definite signs of a development in this direction. Forthcoming publications 
will, therefore, investigate alternative control strategies (including forecasts). 
First results on how to decrease the instability of supply chains can be found 
in Refs. |H[T71[EE]. 
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6.1 Advantages, Extensions, and Potential Applica- 
tions 

As the variables in our dynamical model of supply and production networks 
are operational and measurable, the model can be tested and calibrated with 
empirical data. Moreover, it is flexible and easy to generalize. For this reason, 
it can be adapted to various applications. Our approach can be related to 
microscopic considerations such as queueing theory or event-driven (Monte- 
Carlo) simulations of production processes, but, as it focusses on the average 
dynamics, it is numerically much more efficient and, therefore, suitable for 
on-line control. Nevertheless, the formulas can be extended by noise terms to 
reflect stochastic effects. Our system of coupled differential equations would 
then become a coupled system of stochastic differential equations (Langevin 
equations), where the noise amplitudes would be determined via relationships 
from queueing theory. 

The dynamical theory of supply chains appears to be a promising field 
with many research opportunities. It is not only useful for a deeper under- 
standing of the origin and dynamics of business cycles or for the optimization 
of production processes and supply networks. It could also contribute to the 
further improvement of existing traffic control strategies in urban street net- 
works or to the development of more robust routing algorithms for internet 
traffic. Apart from this, the model can be viewed as a dynamic multi-player 
game, where the individual control functions Wj, reflect the strategies of the 
players b in terms of quantities iVj, which represent the information- feedback 
available to them. Generalizing this idea, the above model of supply net- 
works may serve as a basis for particular kinds of neural networks. We are 
now setting up different projects in these directions, and cooperation is very 
welcome. 
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